Data were loaded from spreadsheets obtained by Utah DWQ:
# tributary data
trib <- read.csv("./Data/Processed/ul_Tribdata_wqp_processed_2020-12-17.csv")
# facility data
fac <- read.csv("./Data/Processed/ul_Facilitydata_wqp_processed.csv")
# watershed spatial data
wshds <- read_sf("./Data/GIS/ulSubWshdsWGS84.shp")
# sampling activity logs
activitylog <- read.csv("./Data/Processed/UL_Trib_SampleActivities_exported 12-17-20.csv")
# groundwater data
gw <- read.csv("./Data/Processed/ul_wqp_gwdata_processed_2021-01-07.csv")
# precipitation and evaporation data
p_et <- read.xlsx("./Data/Processed/ULEFDCFlowBoundaries_WY2006-2018.xlsx",
sheet = "Precip-ET", startRow = 2)
Watersheds were divided into monitored and unmonitored and arranged in clockwise order.
wshds_monitored <- wshds %>%
filter(Monitored_ == "YES")
wshds_monitored$WS_Name <-
factor(wshds_monitored$WS_Name,
levels = c("Lake Mountain - Israel Canyon", "Dry Creek - Saratoga",
"Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough", "Currant Creek"))
levels(wshds_monitored$WS_Name)
## [1] "Lake Mountain - Israel Canyon" "Dry Creek - Saratoga"
## [3] "Lehi Spring Creek" "American Fork River"
## [5] "Timp SSD" "Lindon Drain"
## [7] "Powell Slough Major" "Provo River"
## [9] "Mill Race" "Hobble Creek"
## [11] "Dry Creek - Spanish Fork" "Spanish Fork River"
## [13] "4000 South Drain Spanish Fork" "Benjamin Slough"
## [15] "Currant Creek"
Tributary dataset was filtered for:
Flow data were converted to cubic feet/sec (cfs) when units originally appeared as gallon/min, cubic meters/sec, and million gallons/day.
trib$ActivityStartDate <- as.Date(trib$ActivityStartDate)
# unique(trib$CharacteristicName)
trib$ResultMeasureValue <- as.numeric(as.character(trib$ResultMeasureValue))
## Warning: NAs introduced by coercion
trib_wq <- trib %>%
# Set nondetects (qualifier code U) as 1/2 MDL.
# Qualifier code J is above detection but below reporting. Values retained as-is.
# REVISIT: What does qualifier code Q mean? One sample for ammonium at 0.062
mutate(ResultMeasureValue = case_when(MeasureQualifierCode == "U" ~ MethodDetectionLevel/2,
TRUE ~ ResultMeasureValue)) %>%
filter(CharacteristicName %in%
c("Ammonia-nitrogen as N", "Carbon", "Nitrate as N", "Organic carbon",
"Ammonium", "Carbonate", "Inorganic carbon", "Nitrate", "Nitrogen",
"Phosphorus", "Ammonia and ammonium",
"Inorganic nitrogen (nitrate and nitrite)", "Nitrite", "Orthophosphate",
"Kjeldahl nitrogen", "Nitrogen, mixed forms (NH3), (NH4), organic, (NO2) and (NO3)",
"Organic Nitrogen", "Ammonia-nitrogen", "Phosphate-phosphorus",
"Total Kjeldahl nitrogen")) %>%
filter(ResultMeasure.MeasureUnitCode %in% c("ueq/L", "mg/l", "mg/l as N", "<NA>", "mg/l as P")) %>%
select(MonitoringLocationIdentifier:ResultSampleFractionText,
OrganizationIdentifier:MonitoringLocationTypeName, LatitudeMeasure,
LongitudeMeasure, ResultMeasureValue, ResultMeasure.MeasureUnitCode,
ResultStatusIdentifier) %>%
mutate(Month = month(ActivityStartDate),
Year = year(ActivityStartDate)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
filter(!is.na(ResultMeasureValue)) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultMeasure.MeasureUnitCode) %>%
unite("Variable", CharacteristicName:ResultSampleFractionText)
trib_flow <- trib %>%
filter(CharacteristicName %in% c("Flow", "Stream flow, instantaneous", "Stream flow, mean. daily")) %>%
droplevels() %>%
#filter(ResultStatusIdentifier != "Provisional") %>% # can filter out provisional data
select(MonitoringLocationIdentifier:ResultSampleFractionText,
OrganizationIdentifier:MonitoringLocationTypeName, LatitudeMeasure,
LongitudeMeasure, ResultMeasureValue, ResultMeasure.MeasureUnitCode,
ResultStatusIdentifier) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultSampleFractionText,
-MonitoringLocationTypeName, -CharacteristicName) %>%
group_by(MonitoringLocationIdentifier, ActivityStartDate, OrganizationIdentifier,
OrganizationFormalName, MonitoringLocationName, LatitudeMeasure,
LongitudeMeasure, ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
summarise(ResultMeasureValue = mean(ResultMeasureValue)) %>%
mutate(Month = month(ActivityStartDate),
Year = year(ActivityStartDate)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
pivot_wider(names_from = ResultMeasure.MeasureUnitCode, values_from = ResultMeasureValue) %>%
# USGS data has both cfs and cms, always duplicated. remove cfs.
select(-'ft3/s')
## `summarise()` regrouping output by 'MonitoringLocationIdentifier', 'ActivityStartDate', 'OrganizationIdentifier', 'OrganizationFormalName', 'MonitoringLocationName', 'LatitudeMeasure', 'LongitudeMeasure', 'ResultMeasure.MeasureUnitCode' (override with `.groups` argument)
colnames(trib_flow)[11:14] <- c("m3.sec", "cfs", "mgd", "g.min")
# convert everything to cfs
trib_flow$mgd <- trib_flow$mgd/0.646317
trib_flow$g.min <- trib_flow$g.min*0.002228017
trib_flow$m3.sec <- trib_flow$m3.sec*35.3147
trib_flow_long <- trib_flow %>%
pivot_longer(cols = c("cfs", "mgd", "g.min", "m3.sec"),
names_to = "orig.units", values_to = "Flow.cfs") %>%
drop_na(Flow.cfs)
# conversions seem to match up
ggplot(trib_flow_long, aes(x = ActivityStartDate, y = Flow.cfs, color = orig.units)) +
geom_point(alpha = 0.5) +
scale_y_log10() +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.8)
## Warning: Transformation introduced infinite values in continuous y-axis
A sample type of “11” indicates a site was not sampled. A sample type of “10” indicates there was no flow. The latter should appear in the flow data as a zero; this was verified by inspecting the two datasets.
Locations and dates when sites were not sampled:
## [1] LOWER S FK PROVO R AT GAGING STATIION
## [2] N FK PROVO R AB CNFL / PROVO R AT WILDWOOD
## [3] CURRANT CK AB MONA RES
## [4] Timpanogos Effluent below constructed duck ponds
## [5] Dry Ck Near Utah Lake-WLA
## [6] DRY CK EAST TRIBUTARY
## [7] AMERICAN FK CK 2.5MI S OF AM FK CITY
## [8] Powell Slough WMA North Outfall to Utah Lake
## [9] Powell Slough WMA South Outfall to Utah Lake
## [10] Saratoga Springs at Cedar Valley
## 10 Levels: AMERICAN FK CK 2.5MI S OF AM FK CITY ... Timpanogos Effluent below constructed duck ponds
## [1] "2011-01-26" "2014-02-19" "2012-01-19" "2012-02-15" "2012-03-21"
## [6] "2012-04-10" "2012-06-20" "2012-07-25" "2011-10-27" "2011-11-16"
## [11] "2011-12-06" "2011-12-07" "2017-05-17" "2017-07-10" "2017-12-20"
## [16] "2017-12-22" "2018-01-10" "2018-02-15" "2018-03-13" "2018-04-19"
## [21] "2018-05-17" "2018-11-20" "2019-05-14" "2019-06-13" "2019-06-17"
## [26] "2019-07-08" "2019-08-12" "2019-09-24" "2019-09-23" "2019-10-08"
## [31] "2019-11-04" "2019-12-10"
Tributary dataset was filtered for:
Flow data were converted to cubic feet/sec (cfs) when units originally appeared as gallon/min, cubic meters/sec, and million gallons/day.
fac$datetime <- as.POSIXct(fac$datetime)
#unique(fac$CharacteristicName)
fac_wq <- fac %>%
# Set nondetects (qualifier code U) as 1/2 MDL.
# Qualifier code J is above detection but below reporting. Values retained as-is.
# REVISIT: What does qualifier code Q mean? One sample for ammonium at 0.062
mutate(ResultMeasureValue = case_when(MeasureQualifierCode == "U" ~ MethodDetectionLevel/2,
TRUE ~ ResultMeasureValue)) %>%
filter(CharacteristicName %in%
c("Ammonia-nitrogen as N", "Carbon", "Nitrate as N", "Organic carbon",
"Ammonium", "Carbonate", "Inorganic carbon", "Nitrate", "Nitrogen",
"Phosphorus", "Ammonia and ammonium",
"Inorganic nitrogen (nitrate and nitrite)", "Nitrite", "Orthophosphate",
"Kjeldahl nitrogen", "Nitrogen, mixed forms (NH3), (NH4), organic, (NO2) and (NO3)",
"Organic Nitrogen", "Ammonia-nitrogen", "Phosphate-phosphorus",
"Total Kjeldahl nitrogen")) %>%
filter(ResultMeasure.MeasureUnitCode %in% c("ueq/L", "mg/l", "mg/l as N", "<NA>", "mg/l as P")) %>%
select(MonitoringLocationIdentifier:LongitudeMeasure, ResultMeasureValue,
ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
mutate(Month = month(datetime),
Year = year(datetime)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
filter(!is.na(ResultMeasureValue)) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultMeasure.MeasureUnitCode) %>%
unite("Variable", CharacteristicName:ResultSampleFractionText)
# match date field with trib
colnames(fac_wq)[2] <- "ActivityStartDate"
fac_flow <- fac %>%
filter(CharacteristicName == "Flow") %>%
#c("Height, gage", "Stream flow, instantaneous", "Stream flow, mean. daily")) %>%
droplevels() %>%
#filter(ResultStatusIdentifier != "Provisional") %>% # can filter out provisional data
select(MonitoringLocationIdentifier:LongitudeMeasure, ResultMeasureValue,
ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultSampleFractionText,
-MonitoringLocationTypeName, -CharacteristicName) %>%
group_by(MonitoringLocationIdentifier, datetime, OrganizationIdentifier,
OrganizationFormalName, MonitoringLocationName, LatitudeMeasure,
LongitudeMeasure, ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
summarise(ResultMeasureValue = mean(ResultMeasureValue)) %>%
mutate(Month = month(datetime),
Year = year(datetime)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
pivot_wider(names_from = ResultMeasure.MeasureUnitCode, values_from = ResultMeasureValue)
## `summarise()` regrouping output by 'MonitoringLocationIdentifier', 'datetime', 'OrganizationIdentifier', 'OrganizationFormalName', 'MonitoringLocationName', 'LatitudeMeasure', 'LongitudeMeasure', 'ResultMeasure.MeasureUnitCode' (override with `.groups` argument)
fac_flow$datetime <- as.Date(fac_flow$datetime)
colnames(fac_flow)[2] <- "ActivityStartDate"
colnames(fac_flow)[11:14] <- c("mgd", "cfs", "g.min", "m3.sec")
# convert everything to cfs
fac_flow$mgd <- fac_flow$mgd/0.646317
fac_flow$g.min <- fac_flow$g.min*0.002228017
fac_flow$m3.sec <- fac_flow$m3.sec*35.3147
fac_flow_long <- fac_flow %>%
pivot_longer(cols = c("cfs", "mgd", "g.min", "m3.sec"),
names_to = "orig.units", values_to = "Flow.cfs") %>%
drop_na(Flow.cfs)
Tributary and facility water quality and flow datasets were combined.
USGS_flow <- readNWISdv(siteNumbers = c("10146400", "10150500",
"10153100", "10163000"),
parameterCd = "00060",
startDate = "2010-01-01")
USGS_flow_processed <- USGS_flow %>%
mutate(WS_Name = case_when(site_no == "10146400" ~ "Currant Creek",
site_no == "10150500" ~ "Spanish Fork River",
site_no == "10153100" ~ "Hobble Creek",
site_no == "10163000" ~ "Provo River")) %>%
select(WS_Name, Date, X_00060_00003)
colnames(USGS_flow_processed) <- c("WS_Name", "ActivityStartDate", "Flow.USGS.cfs")
REVISIT
TO DO: Use area vs. elevation relationship to establish average monthly lake area, with variability quantified.
REVISIT: P and ET only available 2005-2018 from EFDC output. We can only use four years of data if we restrict to 2015 onward. Should we use all years?
Sites located downstream of diversion
Sites that may be backwatered at times
Intermittent sites
Effluent-dominated sites
Sites included in ULWQS SAP
trib_fac_sites <- trib_fac_wq_flow %>%
ungroup() %>%
select(MonitoringLocationIdentifier, MonitoringLocationName,
MonitoringLocationTypeName, LatitudeMeasure, LongitudeMeasure) %>%
distinct() %>%
mutate(SiteType = ifelse(MonitoringLocationTypeName %in%
c("Canal Drainage", "Canal Irrigation", "Canal Transport"), "Canal",
ifelse(MonitoringLocationTypeName %in% c("River/Stream", "Stream", "Stream: Ditch"), "Stream",
ifelse(MonitoringLocationTypeName %in% c("Facility Industrial", "Facility Other"), "Facility",
"Storm Sewer"))))
trib_fac_sites_sf <- st_as_sf(trib_fac_sites, coords = c("LongitudeMeasure", "LatitudeMeasure"),
crs = 4326, dim = "XY") #CRS = WGS84
trib_fac_sites_intersect <- trib_fac_sites_sf %>%
st_intersection(wshds) %>%
left_join(trib_fac_wq_flow, .,)
unique(trib_fac_sites_intersect$MonitoringLocationName[is.na(trib_fac_sites_intersect$WS_Name)])
## [1] Provo River Delta Restoration
## [2] JORDAN R AT UTAH L OUTLET U121 XING
## [3] Powell Slough WMA North Outfall to Utah Lake
## [4] Powell Slough WMA South Outfall to Utah Lake
## [5] Dry Ck Near Utah Lake-WLA
## [6] Unnamed flow at 4000 West @ mouth
## [7] American Fork River at mouth
## [8] Timpanogos WWTP at Utah Lake mouth
## [9] Lindon Drain at Utah Lake Inlet
## [10] Beer Creek at Utah Lake mouth
## [11] Hobble Creek at Utah Lake mouth
## [12] Mill Race at mouth 1.25 miles west of I-15
## [13] Dry Creek North Fork @ Confluence
## [14] Dry Creek South Fork @ Confluence
## 153 Levels: AMERICAN FK CK 2.5MI S OF AM FK CITY ...
unique(trib_fac_sites_intersect$MonitoringLocationIdentifier[is.na(trib_fac_sites_intersect$WS_Name)])
## [1] UTAHDWQ_WQX-4917343 UTAHDWQ_WQX-4994790 UTAHDWQ_WQX-4995210
## [4] UTAHDWQ_WQX-4995230 UTAHDWQ_WQX-4996040 WFWQC_UT-4917712
## [7] WFWQC_UT-4994958 WFWQC_UT-4995043 WFWQC_UT-4995075
## [10] WFWQC_UT-4995210 WFWQC_UT-4995467 WFWQC_UT-4996040
## [13] WFWQC_UT-4996096 WFWQC_UT-4996536 WFWQC_UT-4996047
## [16] WFWQC_UT-4996048
## 167 Levels: UTAHDWQ_WQX-4917343 UTAHDWQ_WQX-4993000 ... WFWQC_UT-4996020
unique(trib_fac_sites_intersect$WS_Name)
## [1] NA "Lehi Drain 2"
## [3] "Lake Mountain - Israel Canyon" "Dry Creek - Saratoga"
## [5] "Lehi Spring Creek" "American Fork River"
## [7] "Timp SSD" "Lindon Drain"
## [9] "Powell Slough Major" "Currant Creek"
## [11] "Benjamin Slough" "Spanish Fork River"
## [13] "Dry Creek - Spanish Fork" "Provo Bay 3"
## [15] "Provo Bay 5" "Hobble Creek"
## [17] "Mill Race" "Provo River"
## [19] "Vineyard Drain 6" "4000 South Drain Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Powell Slough WMA North Outfall to Utah Lake"] <- "Powell Slough Major"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Powell Slough WMA South Outfall to Utah Lake"] <- "Powell Slough Major"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Dry Ck Near Utah Lake-WLA"] <- "Dry Creek - Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Mill Race at mouth 1.25 miles west of I-15"] <- "Mill Race"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Hobble Creek at Utah Lake mouth"] <- "Hobble Creek"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Beer Creek at Utah Lake mouth"] <- "Benjamin Slough"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Lindon Drain at Utah Lake Inlet"] <- "Lindon Drain"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Timpanogos WWTP at Utah Lake mouth"] <- "Timp SSD"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"American Fork River at mouth"] <- "American Fork River"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Unnamed flow at 4000 West @ mouth"] <- "4000 South Drain Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"JORDAN R AT UTAH L OUTLET U121 XING"] <- "Jordan River"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"DRY CK EAST TRIBUTARY"] <- "Dry Creek - Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Drainage Canal 0.5 mile bl I-15 at about 2500 West, Springville"] <- "Dry Creek - Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"PROVO STATION 6-WLA"] <- "Mill Race"
trib_fac_sites_wsds <- trib_fac_sites_intersect %>%
ungroup() %>%
select(MonitoringLocationIdentifier, OrganizationIdentifier:LongitudeMeasure,
SiteType:Shape_Area) %>%
distinct() %>%
arrange(WS_Name)
# write.csv(trib_fac_sites_wsds, file = "./Data/Processed/trib_fac_sites_initial.csv",
# row.names = FALSE)
trib_fac_sites_counts <- trib_fac_sites_intersect %>%
ungroup() %>%
group_by(MonitoringLocationIdentifier, OrganizationIdentifier, OrganizationFormalName,
MonitoringLocationName, MonitoringLocationTypeName, LatitudeMeasure,
LongitudeMeasure, SiteType, OBJECTID_1, OBJECTID, TARGET_FID,
Type, wsa_sqkm, WS_Name, Monitored_, WS_NAME_1, WS_NAME_12,
Shape_Leng, Shape_Le_1, Shape_Area) %>%
tally() %>%
ungroup() %>%
select(MonitoringLocationIdentifier, n)
# trib_fac_sites_SD <- read.csv("./Data/Processed/trib_fac_sites.csv")
#
# trib_fac_sites_combined <- full_join(trib_fac_sites_SD, trib_fac_sites_counts)
# write.csv(trib_fac_sites_combined, file = "./Data/Processed/trib_fac_sites_counts.csv",
# row.names = FALSE)
trib_fac_sites_wsds <- trib_fac_sites_wsds %>%
full_join(., trib_fac_sites_counts)
# write.csv(trib_fac_sites_wsds, file = "./Data/Processed/trib_fac_sites_counts.csv",
# row.names = FALSE)
trib_fac_sites_intersect_downstream <- trib_fac_sites_intersect %>%
filter(MonitoringLocationIdentifier %in%
c("UTAHDWQ_WQX-4995312", "UTAHDWQ_WQX-4995310", # Currant Creek, duplicate sites at same lat/lon
"UTAHDWQ_WQX-4995465", "WFWQC_UT-4995467", # Benjamin Slough
"UTAHDWQ_WQX-5919910", "WFWQC_UT-4917712", # 4000 South Drain Spanish Fork
"UTAHDWQ_WQX-4995578", "WFWQC_UT-4995575", # Spanish Fork River
"UTAHDWQ_WQX-4996040", "WFWQC_UT-4996040", # Dry Creek - Spanish Fork
"UTAHDWQ_WQX-4996042", "UTAHDWQ_WQX-4996044", # Dry Creek - Spanish Fork
"WFWQC_UT-4996100", "UTAHDWQ_WQX-4996100", "WFWQC_UT-4996096",
"UTAHDWQ_WQX-4996275", "WFWQC_UT-4996275", # Hobble Creek
"WFWQC_UT-4996536", "WFWQC_UT-4996540", "UTAHDWQ_WQX-4996540", "UTAHDWQ_WQX-4996566", # Mill Race
"WFWQC_UT-4996680", "UTAHDWQ_WQX-4996680", # Provo River
"UTAHDWQ_WQX-4995230", "WFWQC_UT-4995210", "UTAHDWQ_WQX-4995210", # Powell Slough Major
"WFWQC_UT-4995075", "UTAHDWQ_WQX-4995120", # Lindon Drain
"WFWQC_UT-4995043", "UTAHDWQ_WQX-4995041", "UTAHDWQ_WQX-4995038", # Timp SSD
"WFWQC_UT-4994958", "UTAHDWQ_WQX-4994960", # American Fork River
"WFWQC_UT-4994948", "UTAHDWQ_WQX-4994950", #Lehi Spring Creek
"UTAHDWQ_WQX-4994804", # Dry Creek - Saratoga
"UTAHDWQ_WQX-4994792")) %>% #Lake Mountain - Israel Canyon
mutate(Backwater = ifelse(MonitoringLocationIdentifier %in%
c("WFWQC_UT-4917712", "WFWQC_UT-4994958", "UTAHDWQ_WQX-4995465",
"WFWQC_UT-4995467", "WFWQC_UT-4996536", "UTAHDWQ_WQX-4995210",
"UTAHDWQ_WQX-4995230", "WFWQC_UT-4995210", "WFWQC_UT-4995575"),
"Yes", "No")) %>%
left_join(., USGS_flow_processed)
unique(trib_fac_sites_intersect_downstream$WS_Name)
## [1] "Lake Mountain - Israel Canyon" "Dry Creek - Saratoga"
## [3] "Lehi Spring Creek" "American Fork River"
## [5] "Timp SSD" "Lindon Drain"
## [7] "Powell Slough Major" "Currant Creek"
## [9] "Benjamin Slough" "Spanish Fork River"
## [11] "Dry Creek - Spanish Fork" "Hobble Creek"
## [13] "Mill Race" "Provo River"
## [15] "4000 South Drain Spanish Fork"
# check that flow at USGS gage is in line with measured flow at chem sites, when data exist together.
# Some cases in Spanish Fork River where measured flow was near zero but USGS flow was not. Otherwise lines up on 1:1 line.
# Will proceed to replace missing flow with USGS flow in those
# ggplot(subset(trib_fac_sites_intersect_downstream, WS_Name %in%
# c("Hobble Creek", "Spanish Fork River", "Provo River", "Currant Creek")),
# aes(x = Flow.cfs, y = Flow.USGS.cfs, color = WS_Name)) +
# geom_point() +
# geom_abline(slope = 1, intercept = 0)
trib_fac_sites_intersect_downstream$Flow.cfs[is.na(trib_fac_sites_intersect_downstream$Flow.cfs)] <-
trib_fac_sites_intersect_downstream$Flow.USGS.cfs[is.na(trib_fac_sites_intersect_downstream$Flow.cfs)]
# remove zero flow samples from watersheds that do not go dry
# then, calculate daily loads
trib_fac_sites_intersect_downstream <- trib_fac_sites_intersect_downstream %>%
mutate(Flow.cfs = case_when(WS_Name == "Lehi Spring Creek" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Timp SSD" & OrganizationIdentifier == "WFWQC_UT" &
Flow.cfs == 0 ~ NA_real_,
WS_Name == "Powell Slough Major" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Provo River" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Mill Race" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Hobble Creek" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Dry Creek - Spanish Fork" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Spanish Fork River" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "4000 South Drain Spanish Fork" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Benjamin Slough" & Flow.cfs == 0 ~ NA_real_,
TRUE ~ Flow.cfs)) %>%
mutate(load.TN.kgd = TN * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TDN.kgd = TDN * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TP.kgd = TP * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TDP.kgd = TDP * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TOC.kgd = TOC * Flow.cfs * 28.3168 * 86400 / 1000000,
load.DOC.kgd = DOC * Flow.cfs * 28.3168 * 86400 / 1000000)
trib_fac_sites_intersect_downstream$load.TN.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TDN.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TP.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TDP.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TOC.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.DOC.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$SiteType[is.na(trib_fac_sites_intersect_downstream$SiteType)] <- "Not specified"
trib_fac_sites_downstream <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
select(MonitoringLocationIdentifier, MonitoringLocationName, OrganizationIdentifier,
MonitoringLocationTypeName, LatitudeMeasure, LongitudeMeasure, WS_Name) %>%
distinct()
trib_fac_sites_downstream_sf <- st_as_sf(trib_fac_sites_downstream, coords = c("LongitudeMeasure", "LatitudeMeasure"),
crs = 4326, dim = "XY") #CRS = WGS84
trib_fac_sites_downstream_sf$WS_Name <-
factor(trib_fac_sites_downstream_sf$WS_Name,
levels = c("Lake Mountain - Israel Canyon", "Dry Creek - Saratoga",
"Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough", "Currant Creek"))
trib_fac_sites_intersect_downstream$WS_Name <-
factor(trib_fac_sites_intersect_downstream$WS_Name,
levels = c("Lake Mountain - Israel Canyon", "Dry Creek - Saratoga",
"Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough", "Currant Creek"))
# Create two datasets: one starting in 2010, one starting in 2015. The later period of record is the main dataset.
trib_fac_sites_intersect_downstream_early <- trib_fac_sites_intersect_downstream
trib_fac_sites_intersect_downstream <- trib_fac_sites_intersect_downstream %>%
filter(Year >= 2015)
The top graph depicts all tributary and facility monitoring sites. The second graph depicts each monitored watershed individually from the NW side of the lake moving clockwise, with downstream monitoring sites marked. The third graph depicts all watersheds with only downstream monitoring sites marked.
There is only one inflow monitored in Lake Mountain - Israel Canyon. Estimates are based on this site.
Flow
Nutrient loads
Calculate mean monthly loading
There is only one inflow monitored in Dry Creek - Saratoga. Estimates are based on this site.
Flow
Nutrient loads
Calculate mean monthly loading
There are two sites in Lehi Spring Creek. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Calculate mean monthly loading
There are two sites in American Fork River. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Calculate mean monthly loading
Flow
Nutrient loading
Calculate mean monthly loading
REVISIT East side tributary (UTAHDWQ_WQX-4995041) has measurable concentrations but no flow (one near-zero flow, one flow duplicated from 4995038). Scott’s recommendation was to average concentrations at 4995041 and 4995038 and apply them to the flow at 4995038, but this assumes that the volume-weighted average is the same as the arithmetic average (i.e., that the volumes at both sites are the same). I’m not sure if this assumption is valid. For now, I will ignore the east side tributary and calculate loads based on the sites that have flow data.
There are two sites in Lindon Drain. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Calculate mean monthly loading
Flow
Nutrient loading
Calculate mean monthly loads
Powell Slough has both North and South outfalls. Need to add together load from North and South to compute total load. WFWQC samples only at the north outfall.
There are two site codes in Provo River (one DWQ and one WFWQC), but these are in the same location. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Calculate mean monthly loading
Flow
Nutrient loading
Calculate mean monthly loads
4996536 represents the sum total of loading from the watershed. 4996540 and 4996566 should be added together to represent total flow.
Flow
Nutrient loading
Calculate mean monthly loads
4996096 and 4996100 are equivalent sites. 4996275 is a separate flow representing Springville and represents its own load.
Flow
Nutrient loading
Calculate mean monthly loads
4996040 is the primary site in this watershed and represents the sum total of loads. 4996042 and 4996044 represent equivalent sites for the East tributary but do not include the south fork. East tributary sites will be removed and the primary site used for loads.
There are two sites in Spanish Fork River. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Calculate mean monthly loading
There are two sites in 4000 South Drain Spanish Fork. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Calculate mean monthly loading
There are two sites in Benjamin Slough. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Calculate mean monthly loading
Two sites, UTAHDWQ_WQX-4995310 and UTAHDWQ_WQX-4995312, are located at the same coordinates. These sites are pooled in load calculations, as they represent the same location.
Flow
Nutrient loads
Calculate mean monthly loading
loading_AllWatersheds <- rbind(
loading_LakeMountainIsraelCanyon, loading_DryCreekSaratoga,
loading_LehiSpringCreek, loading_AmericanForkRiver, loading_TimpSSD,
loading_LindonDrain, loading_PowellSloughMajor, loading_ProvoRiver,
loading_MillRace, loading_HobbleCreek, loading_DryCreekSpanishFork,
loading_SpanishFork, loading_4000SouthDrain, loading_BenjaminSlough,
loading_CurrantCreek)
# calculate difference between UTAHDWQ and WFWQC load estimates
loading_AllWatersheds_differences <- loading_AllWatersheds %>%
filter(WS_Name %in% c("Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough")) %>%
select(-(load.TN.kgd:load.DOC.kgd)) %>%
group_by(WS_Name, Month) %>%
summarise(load.TN.difference = load.TN.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TN.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TDN.difference = load.TDN.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TDN.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TP.difference = load.TP.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TP.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TDP.difference = load.TDP.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TN.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TOC.difference = load.TOC.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TOC.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.DOC.difference = load.DOC.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.DOC.tonmo[OrganizationIdentifier == "WFWQC_UT"])
## `summarise()` regrouping output by 'WS_Name', 'Month' (override with `.groups` argument)
loading_AllWatersheds_annual <- loading_AllWatersheds %>%
group_by(WS_Name, OrganizationIdentifier) %>%
summarise(Flow.cfs = sum(Flow.cfs),
load.TN.tonyr = sum(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonyr = sum(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonyr = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonyr = sum(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonyr = sum(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonyr = sum(load.DOC.tonmo, na.rm = TRUE))
## `summarise()` regrouping output by 'WS_Name' (override with `.groups` argument)
loading_Total_annual <- loading_AllWatersheds_annual %>%
filter(OrganizationIdentifier == "UTAHDWQ_WQX") %>%
ungroup() %>%
summarise(Flow.cfs = sum(Flow.cfs),
load.TN.tonyr = sum(load.TN.tonyr, na.rm = TRUE),
load.TDN.tonyr = sum(load.TDN.tonyr, na.rm = TRUE),
load.TP.tonyr = sum(load.TP.tonyr, na.rm = TRUE),
load.TDP.tonyr = sum(load.TDP.tonyr, na.rm = TRUE),
load.TOC.tonyr = sum(load.TOC.tonyr, na.rm = TRUE),
load.DOC.tonyr = sum(load.DOC.tonyr, na.rm = TRUE))
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TN.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TN load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 34 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TDN.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TDN load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 105 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TP.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TP load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 24 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TDP.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TDP load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 34 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TOC.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TOC load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 105 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.DOC.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "DOC load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 105 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds,
aes(x = WS_Name, y = load.TN.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
facet_wrap(vars(OrganizationIdentifier)) +
labs(x = "", y = "TN loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 35 rows containing missing values (position_stack).
ggplot(subset(loading_AllWatersheds, OrganizationIdentifier == "UTAHDWQ_WQX"),
aes(x = WS_Name, y = load.TDN.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "TDN loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
ggplot(loading_AllWatersheds,
aes(x = WS_Name, y = load.TP.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
facet_wrap(vars(OrganizationIdentifier)) +
labs(x = "", y = "TP loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 24 rows containing missing values (position_stack).
ggplot(loading_AllWatersheds,
aes(x = WS_Name, y = load.TDP.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
facet_wrap(vars(OrganizationIdentifier)) +
labs(x = "", y = "TDP loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 37 rows containing missing values (position_stack).
ggplot(subset(loading_AllWatersheds, OrganizationIdentifier == "UTAHDWQ_WQX"),
aes(x = WS_Name, y = load.TOC.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "TOC loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
ggplot(subset(loading_AllWatersheds, OrganizationIdentifier == "UTAHDWQ_WQX"),
aes(x = WS_Name, y = load.DOC.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "DOC loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
TO DO: Change y axis to volumes after determining monthly lake area